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Abstract 

Background: According to the growth pattern of bamboo, sympodial bamboo and monopodial bamboo are considered as 
two mainly kinds of bamboo. They have different phenotypes and different characteristics in developmental stage. Much 
attention had been paid on the study of bamboo cultivation, processing, physiology, biochemistry and molecular biology, 
which had made great progresses in the last decade, especially for the highlighted achievement of the bamboo genomics. 
However, there is no information available on concerning comparative profiling of miRNAs between sympodial bamboo 
and monopodial bamboo, which might play important roles in the regulation of bamboo development. 

Methodology/Principal Findings: We identified the profiles of small RNAs using leaf tissues from one sympodial bamboo 
i.e. moso bamboo {Phyllostachys edulis) and another monopodial bamboo i.e. ma bamboo (Dendrocalamus latiflorus). The 
result showed that there were 19,295,759 and 11,513,888 raw sequence reads, in which 92 and 69 conserved miRNAs, as 
well as 95 and 62 novel miRNAs were identified in moso bamboo and ma bamboo, respectively. The ratio of high conserved 
miRNA families in ma bamboo is more than that in moso bamboo. In addition, a total of 49 and 106 potential targets were 
predicted in moso bamboo and ma bamboo, respectively, in which several targets for novel miRNAs are transcription 
factors that play important roles in plant development. More importantly, annotation of differentially expressed target 
genes was performed based on the analysis of pathway and gene ontology terms enrichment. 

Conclusions/Significance: This study provides the first large-scale sight of discovery and comparative characterization of 
miRNAomes between two representative bamboos belonged to sympodial bamboo and monopodial bamboo, respectively. 
Although it will be necessary to validate the function of miRNAs through more experimental research in further, these 
results lay a foundation for unraveling the miRNA-mediated molecular processes in different kinds of bamboo. 



Citation: Zhao H, Wang L, Dong L, Sun H, Gao Z (2014) Discovery and Comparative Profiling of microRNAs in Representative Monopodial Bamboo {Phyllostachys 
edulis) and Sympodial Bamboo {Dendrocalamus latiflorus). PLoS ONE 9(7): e102375. doi:10.1371/journal.pone.0102375 

Editor: Sue Cotterill, St. Georges University of London, United Kingdom 

Received March 10, 2014; Accepted June 17, 2014; Published July 11, 2014 

Copyright: © 2014 Zhao et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits 
unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. 

Funding: This work was financially supported under the introduction project of "948" by the State Forestry Administration of China (No. 201 1-4-55) and the 
fundamental research fund of International Center for Bamboo and Rattan (No. 1632013009). The funders had no role in study design, data collection and 
analysis, decision to publish, or preparation of the manuscript. 

Competing Interests: The authors have declared that no competing interests exist. 

* Email: gaozhimin@icbr.ac.cn 

9 These authors contributed equally to this work. 



Introduction/Background 

MicroRNAs (miRNAs) are endogenous, non-coding small 
RNAs of approximately 22 nucleotides that regulate the flow of 
genetic information by controlling the translation or stability of 
mRNAs. It has been estimated that miRNAs account for ~ 1 % of 
predicted genes in higher eukaryotic genomes, and that up to 10- 
30% of genes may be regulated by miRNAs, which exist widely in 
animals, plants and some viruses. Since the discovery of the first 
miRNA in Caenorhabditis elegans, a large number of miRNAs have 
been identified by experimental cloning or bioinformatics methods 
with the development of biotechnology. More importantly, with 
the increasing development of high-throughput technology, the 
integrated analysis based on miRNA sequences in the evolutionary 
process provides a significant number of evidences that miRNAs 



evolved independently in different strains and generally were 
conserved in the process of evolution [1,2]. However, some data 
and studies associated function with evolution for bamboo 
miRNAs have been rare till now. 

As a tribe of flowering and evergreen perennial monocots, 
classified in the subfamily Bambusoideae within the grass family 
Poaceae that includes rice, maize, wheat and other cereals [3], 
bamboo represents one of the fastest-growing plants and one of the 
most important non-timber forest resources in the world [4,5]. 
More importantly, bamboo is of notable economic and environ- 
mental significance, due to being used for versatile raw products. 
About 2.5 billion people depend on bamboo economically, and 
the international trade in bamboo amounts to over 2.5 billion US 
dollars per year [6]. 
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Bamboo has rhizomes and roots underground that lives 
exclusively in forest and grows into large woody culms. There 
are two mainly patterns of bamboo for their growth, namely as 
sympodial bamboo and monopodial bamboo. Monopodial bam- 
boo grows fast, through their roots or rhizome, which can spread 
widely underground and send up new culms to break through the 
surface. Therefore, this kind was named as running bamboo 
vividly. Monopodial bamboo ranges in all shapes and sizes, from a 
modest 30 cm {i.e. pygmy bamboo) to a mighty 18-20 m in height 
[i.e. moso bamboo). As one of significant bamboo in monopodial 
bamboo, moso bamboo (Phyllostachys edulis) is the most important 
economic one with many advantages such as fast growing rate, 
high yield, extensive use, short crucial period formation and strong 
regeneration capacity. 

Sympodial bamboo grows like a tussock, and their new clums 
emerged around the outside edge of clump, which was apdy 
named as clumping bamboo. Due to unique root and rhizome 
structures, sympodial bamboo produces a tight cluster of culms 
and spreads only a couple centimeters out from the base each 
season. They range in height from less than 20 cm to more than 
30 m, depending on different species of bamboo. Distinguished 
from monopodial bamboo, sympodial bamboo have some special 
features, such as lower spreading rate and u-shaped rhizome 
making new culms next to original plant. As one of the most 
popular and valuable bamboo species in sympodial bamboo, ma 
bamboo (Dendrocalamus latiflorus) fascinated people because of its 
evergreen color and delicious shoot, which was widely cultivated in 
southern China. 

To facilitate omics studies in bamboo and help unveil the 
miRNAs feature of bamboo, small-RNA libraries of leaf tissues 
from moso bamboo and ma bamboo, respectively, were 
constructed and sequenced. In an attempt to identify new 
conserved miRNAs and novel miRNAs in bamboo, a total of 
7,450 known miRNAs from 70 plants species from miRBase 
served as a reference data set. Moreover, the potential target genes 
of novel miRNAs were predicted based on the high-throughput 
technology and bioinformatics analysis. The profiling of micro- 
RNAs from moso bamboo and ma bamboo were compared. 
Ultimately, these results will play the necessary role in under- 
standing the regulation of miRNAs in the development of 
monopodial bamboo and sympodial bamboo. 

Materials and Methods 

Plant material and RNA isolation 

Divided seedlings of moso bamboo and ma bamboo were potted 
in our laboratory under a regime of 1 6 h light and 8 h darkness at 
25°C, with a light intensity of 200 umol-m -2 -s -1 and a relative 
humidity of 75%. The leaf tissues were used as materials for this 
study. We chose the third piece of new functional leaf (blade tissue 
only) from the top of the branch, which could be considered as a 
juvenile leaf. The leaves were collected and quickly frozen in liquid 
nitrogen. Total RNA was isolated from leaf tissues using the Trizol 
reagent (Invitrogen, Carlsbad, CA, and USA), according to the 
manufacturer's instructions strictly. 

Small RNA high-throughput sequencing 

The small RNA library construction and Solexa sequencing 
were carried out at BGI-Shenzhen (Shenzhen, China) using the 
standard small RNA sample preparation protocol (Illumia) as 
described by Hafner et al. [7]. Briefly, small RNAs of 15-30 nt in 
length were first isolated from the total RNA through 15% TBE 
urea denaturing polyacrylamide gels, which were ligated to the 
small RNA-5' adaptor, and then 3' adaptor was ligated to the 



small RNA-5' adaptor, followed by reverse transcription into 
cDNAs. These cDNAs were amplified by PCR and subjected to 
Solexa sequencing. 

Sequence data analysis 

After removing low quality reads and trimming adapter 
sequences, initial reads based on Solexa sequencing were 
processed by summarizing data production, evaluating sequencing 
quality, calculating the length distribution of small RNA reads. 
Moreover, small RNAs ranging from 18-30 nt were collected and 
used for further analyses. The raw reads are available in the NCBI 
SRA database under the accession number SRX480448. For 
analyzing conserved and novel miRNA, firstly, the clean reads, 
which mapped to exons or other non-coding RNAs such as rRNA, 
snRNA, snoRNA and tRNA, are filtered by BLAST against the 
coding sequences of moso bamboo genome, the Rfam database 
(http://rfam.sanger.ac.uk/) [8] and the GenBank non-coding 
RNA database (http://www.ncbi.nlm.nih.gov/) to discard coding 
gene, non-coding RNA sequences. Subsequently, unique sRNAs 
were aligned with plant mature miRNAs in miRBase Release 20. 
After rigorous screening, sRNA sequences without more than 
three mismatched bases were selected by BLAST searching against 
miRBase. Then, the remaining reads were used to map the 
genome of moso bamboo. Sequences with a tolerance of two 
mismatches were retained for miRNA prediction. RNAfold 
(http://www.tbi.univie.ac.at/RNA/) [9] was used for secondary 
structure prediction (hairpin prediction) of individual mapped 
miRNAs, using the default folding conditions to identify conserved 
miRNAs in bamboo. Lastly, sequences without defined as the 
conserved miRNAs were termed as novel miRNAs. 

RPKM (Reads Per Kilobase per Million) is a method of 
quantifying gene expression from RNA sequencing data by 
normalizing for total read length and the number of sequencing 
reads [10]. Based on the same principle, RPKM is suitable the 
analysis of miRNA as well. The equation for calculating RPKM is 
in the following. 



- C: the number of one reads mapped genome. 

- N: the number of total reads mapped genome. 

- L: the length of reads mapped genome. 

In addition, to reasonably compare the miRNA abundance in 
libraries of moso bamboo and ma bamboo, the value of RPKM 
were calculated to seek the normalization of the difference 
between two samples. As a unique exact match, the information 
that a read mapped only once in the genome with no mismatch 
was also provided. 

Potential target sequences for the novel miRNAs were predicted 
on the basis of the psRNATarget program (http:/ /plantgrn.noble. 
org/psRNATarget/) with default parameters [11]. As the 
requirement of program, the newly identified miRNA sequences 
were used as custom miRNA sequences, while coding sequences in 
moso bamboo genome were used as custom plant databases. All 
predicted target genes were evaluated by the scoring system and 
criteria defined in a previous report [12]. 

To obtain custom enrichments, Gene Ontology (GO) terms 
enrichment was calculated by Ontologizer [13], which is a java 
application that can be used to perform statistical analysis for 
enrichment of GO terms in set of genes or proteins based on the 
one-sided Fisher's exact test, the novel parent-child method, and 
topology-based algorithms. As the inputs for Ontologizer, the four 
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Figure 1. Length distribution of bamboo small RNA. 
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files were needed, including one GO file with the format of open 
biomedical ontologies, one list of all gene IDs of moso bamboo, 
one list of gene IDs that will be analyzed and one association file of 
moso bamboo. As parameter settings, term for term analysis with 
model-based gene set was used. GO terms analysis was performed 
separately on the set of target genes in sample of moso bamboo 
and ma bamboo. Besides, GO categories with p-values less than 
0.05 were considered significant. Lastly, the biological interpreta- 
tion of the target genes of differential miRNAs was completed 
using KEGG (Kyoto Encyclopedia of Genes and Genomes) 
pathway analysis. 

Expression analysis of miRNAs by qRT-PCR 

To detect the expression level of miRNA in bamboo, stem-loop 
qRT-PCR was employed in this study. Primers were designed 
specifically based on individual miRNA using previous method 
[14]. U6 snRNA was selected as internal control [15]. More 
detailed information of primers was supported in File SI. 

On the basis of leaf samples in moso bamboo and ma bamboo, 
cDNAs were synthesized from total RNA with the miRNA-specific 
stem-loop RT primer according to reference [16]. Subsequently, 
qRT-PCR was carried out using a SYBR Green I Master Kit 
(Roche, Germany) on a QTOWER2.2 Real-Time PCR System 
(Analytik Jena). The final volume was 10 ul, containing 5.0 jllI 
2 xSYBR Premix Ex Taq, 0.2 ul of each primer (10 uM), 0.8 ul of 
cDNA and 3.8 ul of nuclease-free water. All reactions were 
repeated three times and the amplification was conducted as 
follows: initial denaturation at 95°C for 10 min, followed by 50 
cycles at 95°C for 10 s, and 62°C for 10 s. There was three 
biological experiments. Expression levels were normalized to that 
of the internal control, and the relative value was calculated using 
previous method [17]. 

MiRNA ioslation by PCR and sequencing 

To prove the existence of miRNAs in bamboo, stem-loop 
primers were also directly used to amplify the mature miRNA 
sequences using specific templates from moso bamboo and ma 
bamboo respectively. Moreover, the primers designed for precur- 
sors (File S2) were used to amplify with DNA templates of moso 
bamboo and ma bamboo respectively. All the amplified products 
were sequenced and analyzed. 



Results and Discussion 

Summary of small RNA library sequencing 

Based on Solexa sequencing of small RNA libraries from moso 
bamboo and ma bamboo (File S3), 19,295,759 and 11,513,888 
raw reads were produced, respectively. After discarding low 
quality sequences and collapsed raw data, there were 8,089,795 
(moso bamboo) and 10,593,305 (ma bamboo) clean reads as well 
as 2,653,717 (moso bamboo) and 6,320,379 (ma bamboo) unique 
sequences. Subsequendy, the sequences, which mapped to coding 
sequences of moso bamboo genome and non-coding RNAs such 
as tRNAs, rRNAs, siRNAs, snRNA, and snoRNA, were removed. 
The statistics of small RNA library sequencing were shown in 
Table 1. As an important feature of the size profile, the 
distribution of small RNA length was summarized in Figure 1. 
The general distribution profile of moso bamboo was similar to 
that of ma bamboo. The most of small RNAs were ranged from 
20 nt to 24 nt in length, which were highly consistent with those of 
small RNAs for known function [18]. As a maximum peak, small 
RNA with 24 nt in length was also highly agreed with those small 
RNAs in other Poaceae plants based on Solexa sequencing 
technology [19-21]. 

Identification of conserved miRNAs 

We identified 92 and 69 conserved miRNAs, belonged to 53 
families and 35 families for moso bamboo and ma bamboo, 
respectively, based on sequence alignment by BLAST against the 
latest miRBase v20, in which miRNAs of bamboo were no 
available, and a series of strict filtering criteria (File S4). 

Normalized analysis of miRNA sequence between moso 
bamboo and ma bamboo were employed. The number of mapped 
reads, the length of mapped reads and the number of total mapped 
reads were considered to make it possible to compare RNA 
expression levels from different samples. In addition, the 
previously reports shown miRNAs with high sequencing frequen- 
cies played some key regulatory functions in maintaining biological 
process [22]. Consequently, the read counts for known miRNA 
families were evaluated in two samples. The result indicated small 
RNA families of moso bamboo were mapped to numerous 
conserved miRNA families, while those of ma bamboo were only 
concentrated in some conserved miRNA families. The reads with 
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RPKM>5000 in two samples were focused on the abundance of 
conserved miRNA families, which accounted for a high proportion 
(~99%) of the total conserved miRNA both in ma bamboo and 
moso bamboo. Apart from high abundant miRNA families, one 
must bear in mind that, there are a large number of low abundant 
miRNA families, indicating that these miRNAs were expressed at 
a low level. Another feature is that some conserved miRNAs were 
identified only in one sample. For example, miR396, miR397, 
miR1432 and miR7748 were identified only in moso bamboo, 
wherever miRl 70 was only in ma bamboo. This may be explained 
by the fact that the two kinds of bamboo had their own diversely 
biological and developmental features in the long progress of 
evolution. However, the previously study on the miRNA of ma 
bamboo demonstrated that miR396 family was identified [12]. It 
was may be caused by tissue and process specific miRNA in two 
different samples of ma bamboo [23,24]. Taken together, as 
specific conserved miRNAs in moso bamboo, miR397, miR1432 
and miR7748 were focused on using experimental validations 
further. 

The number of members in different conserved miRNA families 
was also analyzed. The number of conserved miRNA families with 
more than 3 members was seven both in moso bamboo and ma 
bamboo. In addition, those with 3 members or less in both of moso 
bamboo and ma bamboo were accounted more than half of the 
total miRNA families. 
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Evolution of miRNA in bamboo 

Due to the important role of miRNA in regulation of plant gene 
expression, researchers have focused on the prediction, identifica- 
tion and functional analysis of miRNA in bamboo. However, the 
evolution of miRNA in bamboo has not been reported. To 
investigate the evolutionary roles of the conserved miRNAs, 
comparisons against conserved miRNAs in plants described by 
Cuperus et al [25] were carried out. The outcome showed almost 
of highly conserved miRNA families of bamboo were found (File 
S5), which further proved the previous arguments that miRNAs 
were highly conserved in sequences [26]. In addition, the previous 
study indicated that miRNA families from miRl 56 to miR408 
were identified as high conserved miRNAs in the plant kingdom 
[27,28]. These miRNAs play essential and conserved functions in 
plant development, such as flower and leaf development. MiRNAs 
families behind miR408 are low conserved or non-conserved, 
which may play roles in more species-specific characteristics in 
plant growth and development. These miRNA families were 
considered as the new generation in evolutional process. As shown 
in Figure 2, comparing conservation for miRNA families between 
two samples, a total of RPKM ranged from miRl 56 to miR408 
accounted for 69.1%, 81.0% in moso bamboo and ma bamboo, 
respectively. Moreover, there is higher ratio of low conserved 
miRNAs families in moso bamboo (3 1 %) than that of ma bamboo 
(19%). It was indicated that the implementation in regulatory 
processes, which played more essential roles, might need more 
involvement of miRNA in ma bamboo leaf than those in moso 
bamboo. 

The previously study demonstrated that miRNA sequences 
have to complement with their target gene transcripts to carry out 
their functions, which suggested that miRNA genes have to co- 
evolve with their target genes [29] . MiRNA genes would amplify 
through duplication events similar to those that drive the evolution 
of protein-coding genes: tandem gene duplications, segmental 
duplications, and chromosomal duplications or polyploidization 
[30,31]. In this study, the different of proportion in low conserved 
miRNAs figures out more miRNAs with low conserved will be 
need to regulate newly emerging genes, which help moso bamboo 
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Figure 2. The ratio of high conserved miRNA between moso bamboo and ma bamboo. 
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adapt to the changing environments better than ma bamboo does 
[32,33]. The rhizome of monopodial bamboo can spread laterally 
grown in soil, and also be differentiated into bamboo away from 
the mother plant, while sympodial bamboo grows in clusters 
within a relatively small range. Therefore, monopodial bamboo 
requires more self-regulation to adapt to the environment than 
sympodial bamboo, which might be consistent with the generation 
of lower conserved miRNAs families. Studies have shown that the 
expression of miRNAs is closely related with the developmental 
stages of plant. Based on the analysis of the previous study on the 
miRNAs of moso bamboo [34], the analysis based on percentages 
of conserved miRNA families in different development of 
internode demonstrated that highly conserved miRNAs may play 
more essential roles in earlier development of moso bamboo due to 
some stress response, while non-conserved miRNAs with low 
expression occupied a dominant position in maintaining the 
regulatory function during mature stage of moso bamboo 



(Figure 3). Moreover, the results based on percentages of 
conserved miRNA families in different positions of same 
development stage illustrated that highly conserved miRNAs from 
different positions of same tissue were possibly stable in one stage. 
To further understand the feature of bamboo, the related studies 
on other species of bamboo will be properly carried out. 

Identification of novel miRNAs 

According to the analysis, we predicted 95 and 62 candidates of 
novel miRNAs with more than 3 reads in moso bamboo and ma 
bamboo, respectively. The mature sequences of novel miRNA 
with RPKM>300 in two samples were shown in Table 2. The 
comprehensive information on novel mature miRNAs was 
provided in File S6. 

The number of novel miRNAs and unique mapped reads of ma 
bamboo were less than those of moso bamboo, probably because 
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Table 3. The target genes of Novel miRNAs with RPKM>300 in moso bamboo and ma bamboo. 



miRNA name 


Target ID 


Inhibition 


Target annotation 


phe-miRC5-1 


PH01000594G0440 


Cleavage 


SPL14 - SBP-box gene family member 


phe-miRC5- 1 


PH01000770G0270 


Cleavage 


SPL17 - SBP-box gene family member 


phe-miRC5- 1 


PH01000150G0460 


Cleavage 


SPL14 - SBP-box gene family member 


phe-miRC5- 1 


PH01000050G0170 


Cleavage 


SPL17 - SBP-box gene family member 


phe-miRC5-l 


PH01003178G0220 


Cleavage 


SPL2 - SBP-box gene family member 


phe-miRC5- 1 


PH01000002G1660 


Cleavage 


SPL2 - SBP-box gene family member 


phe-miRC5-l 


PH01000969G0180 


Cleavage 


SPL16 - SBP-box gene family member 


phe-miRC5- 1 


PH01002673G0070 


Cleavage 


SPL16 - SBP-box gene family member 


phe-miRC5-1 


PH01 0001 17G1 390 


Cleavage 


SPL16 - SBP-box gene family member 


phe-miRC5- 1 


PH01003773G0220 


Cleavage 


SPL16 - SBP-box gene family member 


phe-miRC7 


PH01000316G0850 


Cleavage 


PPR repeat containing protein 


phe-miRC9 


PH01000316G0850 


Cleavage 


PPR repeat containing protein 


phe-miRCI 1-1 


PH01000022G2070 


Cleavage 


SAC3/GANP family protein 


phe-miRC10-l 


PH01000316G0850 


Cleavage 


PPR repeat containing protein 


dla-miRCl-1 


PH01000594G0440 


Cleavage 


OsSPL14 - SBP-box gene family member 


dla-miRCl-1 


PH01000770G0270 


Cleavage 


OsSPL17 - SBP-box gene family member 


dla-miRCl-1 


PH01000150G0460 


Cleavage 


OsSPL14 - SBP-box gene family member 


dla-miRCl-1 


PH01000050G0170 


Cleavage 


OsSPL17 - SBP-box gene family member 


dla-miRCl-1 


PH01003178G0220 


Cleavage 


0sSPL2 - SBP-box gene family member 


dla-miRCl-1 


PH01000002G1660 


Cleavage 


OsSPL2 - SBP-box gene family member 


dla-miRCl-1 


PH01000969G0180 


Cleavage 


OsSPL16 - SBP-box gene family member 


dla-miRCl-1 


PH01002673G0070 


Cleavage 


OsSPL16 - SBP-box gene family member 


dla-miRCl-1 


PH01 0001 17G1 390 


Cleavage 


OsSPL16 - SBP-box gene family member 


dla-miRCl-1 


PH01003773G0220 


Cleavage 


OsSPL16 - SBP-box gene family member 


dla-miRC2 


PH01002521G0110 


Cleavage 


acyl-desaturase, chloroplast precursor 


dla-miRC4 


PH01 000371 G0630 


Cleavage 


expressed_protein 


dla-miRC8 


PH01002105G0330 


Cleavage 


OsGrx_C3 - glutaredoxin subgroup 1 


dla-miRC8 


PH01001437G0080 


Cleavage 


regulator of nonsense transcripts 1 


dla-miRC8 


PH01000553G0260 


Cleavage 


regulator of nonsense transcripts 1 


dla-miRC9 


PH01 000371 G0630 


Cleavage 


expressed_protein 


dla-mlRCW 


PH01 000061 G0870 


Cleavage 


tRNA synthetase class 1 


dla-miRCIO 


PH01000087G1230 


Cleavage 


expressed_protein 


dla-miRC6 


PH01 000421 G0420 


Cleavage 


AMP deaminase 


dla-miRC6 


PH01002076G0320 


Cleavage 


expressed_protein 


dla-miRC6 


PH01003497G0060 


Cleavage 


heavy metal-associated domain containing protein 


dla-miRC6 


PH01004252G0100 


Translation 


wall-associated receptor kinase 3 precursor 


dla-miRC6 


PH01003898G0120 


Translation 


zinc finger, ZZ type family protein 


dla-miRC6 


PH01002588G0080 


Cleavage 


expressed_protein 


dla-miRC6 


PH01004366G0020 


Cleavage 


DDT domain containing protein 


dla-miRC6 


PH01000778G0060 


Cleavage 


DDT domain containing protein 


dla-miRC6 


PH01001900G0280 


Cleavage 


expressed_protein 


dla-miRC6 


PH01001125G0220 


Cleavage 


OsSCP19 - Putative Serine Carboxypeptidase homologue 


dla-miRC6 


PH01002299G0250 


Cleavage 


coatomer subunit delta 


dla-miRC7-l 


PH01000594G0440 


Cleavage 


OsSPL14 - SBP-box gene family member 


dla-miRC7-l 


PH01000770G0270 


Cleavage 


OsSPL17 - SBP-box gene family member 


dla-miRC7-l 


PH01000150G0460 


Cleavage 


OsSPL14 - SBP-box gene family member 


dla-miRC7-l 


PH01000050G0170 


Cleavage 


OsSPL17 - SBP-box gene family member 


dla-miRC7-l 


PH01003178G0220 


Cleavage 


OsSPL2 - SBP-box gene family member 


dla-miRC7-l 


PH01000002G1660 


Cleavage 


0sSPL2 - SBP-box gene family member 


dla-miRC7-l 


PH01134604G0010 


Cleavage 


transferase family protein 
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Table 3. Cont. 





miRNA name 


Target ID 


Inhibition 


Target annotation 


dla-miRC7-1 


PH01000197G0200 


Cleavage 


sarcoma antigen NY-SAR-91 


dla-miRC7-1 


PH01001526G0360 


Cleavage 


sister chromatid cohesion 2 


dla-miRC13 


PH01000733G0330 


Cleavage 


HAD-superfamily hydrolase, subfamily IA, variant 3 containing 
protein 


dla-miRC13 


PH01000707G0040 


Cleavage 


oxidoreductase, short chain dehydrogenase/reductase family 
protein 


dla-miRCU 


PH01000031G1370 


Cleavage 


tetratricopeptide repeat domain containing protein 



doi:1 0.1 371 /journal.pone.01 02375.t003 



the genome of moso bamboo was analyzed as required genome. 
Although both of moso bamboo and ma bamboo belonged to 
Bambusoideae, phylogenetic analyses between these had few study 
for the lack of the genome information of ma bamboo. Therefore, 
using the genome of moso bamboo for analyzing miRNA of ma 
bamboo would inevitably produce some biases. To address the 



gap, the whole genome sequencing of ma bamboo need be 
implemented as soon as possible. 

Among the novel miRNAs, /j/ze-miRC 1 - 1 with 106143.04 
RPKM, more than 3-folds of the first-highest in ma bamboo, 
had the highest expression level in moso bamboo. On the basis of 
their frequencies and sequences in the small RNA libraries, 
although the expression levels of these candidates ranged from 



Table 4. List of the KEGG pathways of miRNA target genes affiliated. 





KEGG Pathway 






nu. ui laiytri 

genes 


Target genes 


miRNAs 


Cellular Processes 


Cell growth and death 


Cell cycle 


4 


PH01 00 1 634G02 1 0, PH01 00042 1 G0420, 
PH01003083G0110, PH01000303G0070 


phe-miRC20, 
dla-miRC6, 
dla-miRC34-l, 
dla-miRC34-l 


Cellular Processes 


Cell growth and death 


Oocyte meiosis 


4 


PH01 00 1 634G02 1 0, PH01 00042 1 G0420, 
PH01003083G0110, PH01000303G0070 


phe-miRC20, 
dla-miRC6, 
dla-miRC34-l, 
dla-miRC34-l 


Cellular Processes 


Transport and catabolism 


Peroxisome 


2 


PH01002603G0080, PH01002603G0080 


dla-miRC22, 
dla-miRC35 


Environmental Information 
Processing 


Signal transduction 


TGF-beta signaling 
pathway 


2 


PH01003083G0110, PH01000303G0070 


dla-miRC34-1, 
dla-miRC34-l 


Environmental Information 
Processing 


Signal transduction 


Wnt signaling pathway 


2 


PH01003083G0110, PH01000303G0070 


dla-miRC34-1, 
dla-miRC34-1 


Genetic Information 
Processing 


Folding, sorting and 
degradation 


Protein processing in 
endoplasmic reticulum 


4 


PH01 000063G1 430, PH01 002461 G01 1 0, 
PH01003083G0110, PH01000303G0070 


phe-miRC39, 
dla-miRC28, 
dla-miRC34-l, 
dla-miRC34-l 


Genetic Information 
Processing 


Folding, sorting and 
degradation 


Ubiquitin mediated 
proteolysis 


3 


PH01 000421 G0420, PH01 003083G01 1 0, 
PH01000303G0070 


dla-miRC6, 

dla-miRC34-l, 

dla-miRC34-1 


Genetic Information 
Processing 


Translation 


mRNA surveillance 
pathway 


2 


PH01000553G0260, PH01000553G0260 


dla-miRCB, 
dla-miRC37 


Genetic Information 
Processing 


Translation 


RNA transport 


2 


PH01000553G0260, PH01000553G0260 


dla-miRCa, 
dla-miRC37 


Metabolism 


Amino acid metabolism 


Cysteine and methionine 
metabolism 


1 


PH01001446G0280 


dla-miRC40 


Metabolism 


Amino acid metabolism 


Valine, leucine and 
isoleucine degradation 


1 


PH01001555G0230 


dla-miRC47 


Metabolism 


Carbohydrate metabolism 


Pyruvate metabolism 


1 


PH01002248G0250 


phe-miRC65 


Metabolism 


Energy metabolism 


Carbon fixation in 
photosynthetic organisms 


1 


PH01002248G0250 


phe-miRC65 


Metabolism 


Metabolism of cofactors 
and vitamins 


Ubiquinone and other 

terpenoid-quinone 

biosynthesis 


2 


PH01001360G0010, PH01001360G0010 


dla-miRC21, 
dla-miRC38 



doi:1 0.1 371 /joumal.pone.01 02375.t004 
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a. GO terms enrichment analysis for target genes of moso bamboo 

single-organism process. biological_piocess; 0044699 ^^^^^H^^^^H^^^^H^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^H 

localization; biological_process; GO 0051179 ^^^^^H^^^^^^^^^^^H^^^^h^^^h^^h^^h^^^^^h^^h^^h^^h^hh 

establishment of localization; biological_process; GO 005123-1 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^HI^^BHHMI^^^^^h^h^^^^h^H 

transmembrane transport; biologicaljprooess; GO: 0055085 t^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^M 

tryptopharryl-fRNA aminoacylation; biological_process*, GO 0006436 ^^^^^^^^^^^^h^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^h 

nncroboclv eellularcoiiiponent. GO 004443S ■■^^^^^^^^^^^^^^^^■^^^^^^^^^^^^^^^^^^^^^^^^^^^H 

spindle; cellulai ^component: GO 0005819 ^^^^^^^^^^^^^^^^^^^^^H^^^^^^^^^^^^^^^^^^^^^^^^H 

tryptophan-tRNA ligase activity; in oleculai function; GO 0004830 ■HHHHHHHIHBHHBHBHHaHHHHH 

biotin carboxylase activity; molecular_function. GO 0004075 ^^^^^^■^^^^^^^^^^h^^^^^^^h 

purine ribonucleoside triphosphate binding; molecular_function; GO 0035639 ^^^^^^^^^^^^^^^^^^^^^^^^^h^^^^^^^^^^^^h 

ligase activity, forming carbon-carbon bonds. molecular J:\inetion. GO 001 6885 ^^■^^^^^■^^^^^■^^■^^■^^■^^■^^^^^■^^■B 

single-organism transport; luological_process. GO 0044~65 nnmiiiiiiiiiiiiiiiMiiiiiiiiiiiiiiiiiiiiiiiiiiiiiiM 

binding; arolecular functioir, GO 00054SS BHHHHlllllllllllllllllllllllllllllllllllM 

DNA metabolic process. biological_process. GO 0006259 HHHHHHllllllll^illlllllllllllllllllM 

membrane-bounded organelle, celhilar component; GO 004322" i^is^is^is^is^is^is^is^is^is^is^is^B 

mtra cellulai membrane-bounded organelle, cellular component. GO 0043231 iBHHHHHHHlllllllllllllllllllllllim 

cellular process. biolpgical_process; GO:0009987 iBHHHHHHliMiiiMiiiMiiiMiiH 

DNA polymerase complex. cellular_eoinponent. GO 00425"5 BHIUM 

celhilar_component; celhilar_component. GO 0005575 B1M 

DNA binding; molecuiar_function. GO 00036"" I 
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P-value 



b. GO terms enrichment analysis for target genes of ma bamboo 

methvlauoii: biotogical_process. GO 0032259 
cellular localization. bioIogical_process; GO:0051641 
nucleobase-coiitauiing compound metabolic process. biologieal_proee.ss. GO 0006139 
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AP-rype membrane coat adaptor complex. cellularcomponent. GO 00301 19 

DNA binding; molectilai Junction. GO 0003677 ^^^^^^^^^^^^^h^^^^^^^^^^^^h 
oxidoreducta.se activity, acting on I.2H-OH group of donors. molecii!ai_tunction. GO 001661 4 ^^^^^^^^^^^^^^^^^^z^^^^^^^^^^^^^^^^^^^m 
DNA-methylaaiLif erase activity; molccularrunction, GO: 0009008 ^^^^^^^^^^^^^^^^^^^^^^^^^^ 
cellulai macromolecule catabokc process, biologicaljirocess. GO 0044265 
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nucleobase-contanrrng compound catabolic process; biological_process; GO 0034655 

cellulai catabolic process. bioIogical_proeess. G< I I24S t^^^^^^^^^^^^^^^^^^^^ 

DNA modification; biologicaI_process. GO 0006304 
organic substance catabolic process. biological_process, GO 1901575 ^^^^^^^^^^^^^^^^^^^^^ 
organic cyclic compound binding, molecorlar function; GO 009~159 
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aromatic compound catabolic process, biological_process. GO - 001 9439 
cellulai nitrogen compound catabolic process. bioIogical_proeess. GO ti( i442"n 

cellular nitrogen compound metabolic process. biologieal_piocess. GO 0034641 ^^^^^^^^^^^^^^^ 

cellulai' aromatic compound metabolic process. bioIogieaI_process. GO 0006"25 
catabolic process. biological_process; GO 0009056 
lieterocvcle metabolic process, biological_proeess, GO:00464S3 
nucleic acid binding; molecular fimction; GO 00036~6 
organic cyclic compound catabolic process. biologieal_proeess. GO: 1901361 ^^^^^ m ^^^^ m ^^^ m 
mRNA metabolic process. biological_prooess; v_tO 00160"! ^^mmmmmmi^m^^ 
inRNA catabolic process. biologicaI_process. GO 0006402 
nMmbrane-lxiuuded organelle, cellulai component: GO <XM322~ ■ 
intracellular membrane- bounded organelle, eellularcoiiiponent . GO 0043231 ■ 
RNA catabolic process; biologica process. GO 0006401 ■ 
celluku component . cellular component. GO 00055~5 
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Figure 4. GO terms enrichment analysis for target genes between moso bamboo and ma bamboo. 

doi:10.1371/journal.pone.0102375.g004 

thousands of RPKM to lower, in general, novel miRNA Prediction of miRNA targets 

candidates showed lower expression level compared with most of The acquisition of knowledge for miRNA target genes 

the conserved families. The low abundance of novel miRNAs contributed to insight into the range of miRNA regulation and 
demonstrated that these miRNAs might play a specific and unveiling detailed descriptions of miRNA-target interactions. In 
essential role in certain tissues or developmental stages, and were addition, as an essential step to identify miRNA targets, the 
considered as young miRNAs in terms of evolution [25]. methods based on bioinformatics were widely used [11]. Herein, 

on the basis of the criteria described in the material and methods, 
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Figure 5. Expression analysis of miRNAs in leaves of moso bamboo using qRT-PCR. Error bars representing the standard deviation were 

derived from the three experiments in triplicate. 

doi:10.1371/journal.pone.0102375.g005 



a total of 49 and 1 06 target genes were identified in moso bamboo 
and ma bamboo, respectively. The target genes of novel miRNAs 
with RPKM>300 in two samples were shown in Table 3. The 
comprehensive information on target genes was provided in File 
S7. 

The target genes of the novel miRNAs from bamboo mainly 
concentrated on transcription factor, signal regulation and 
functional protein, among which were a considerable number of 
SPL genes, indicating that these miRNAs might be involved in 
regulating the expression of SPL genes. SPL genes family is a group 
of structurally diverse genes encoding putative transcription factors 
found in photosynthetic organisms [35,36]. The distinguishing 
feature of the SPL gene family is the SBP-box encoding a 
conserved protein domain of 76 amino acids in length. Recent 
studies indicated that SBP-box genes played central roles in 
controlling certain homeostatic processes, as well as in overcoming 
stresses [37,38]. In the previous study of Arabidopsis, several 
targeted SPL genes act redundantly in controlling the juvenile-to- 
adult phase transition, which renders plants sensitive to photope- 
riodic induction of flowering [39,40]. Therefore, the miRNAs 
targeted SPL genes might play important roles in bamboo leaves, 
which response to environmental factors especially for light 
induction. 

Annotation of differentially expressed target genes 

Functional annotations were carried out to investigate which 
processes are differentially regulated by miRNA in moso bamboo 
and ma bamboo. To infer significant associations between the sets 
of targets predicted by individual miRNAs and specific biochem- 
ical pathways, a total of 20 different metabolic pathways were 
found, including 5 metabolism, 8 cellular processes, 7 genetic 
information processing (Table 4), which indicated that miRNA 
played essential roles in various biological process. 

Besides pathway analysis, another functional annotation was 
performed via the analysis of GO, which is the de facto standard in 
gene functionality description and is used widely in functional 
annotation and enrichment analysis [41]. Genes from a large 



amount of organisms have been annotated to GO terms. The 
widespread applications based on GO terms enrichment analysis 
were the identification of annotation-enriched GO terms in a list 
of genes with some similar characteristics in biology. These terms 
were often considered as some representations with the outstand- 
ing biological features of the genes in the study of set [13]. We 
performed GO analysis on the result of prediction of miRNA 
targets using model-based gene set analysis (MGSA), which 
analyzes all GO terms at once by embedding them in a Bayesian 
network in order to provide high-level, summarized views of core 
biological processes. 

In order to identify significant target genes, the GO term for 
target genes with p-value <0.05 were selected and shown in 
Figure 4. GO term enrichment analysis identified a total of 20 and 
32 GO terms with high significant in moso bamboo and ma 
bamboo, respectively. The most highly ranked term of the target 
genes with higher significant in moso bamboo was DNA binding 
(GO: 0003677), which included 17 GO terms from the total 100 
GO terms in the set of the moso bamboo. On the other hand, the 
most highly ranked term for the genes with significant in ma 
bamboo was cellular component (GO: 0005575) including 19 ones 
of the total 69 GO terms in the set of ma bamboo. Moreover, 
there were 4 GO terms overlapped (GO: 0003677, GO: 0005575, 
GO: 0043231 and GO: 0043227) in the total of significantly GO 
terms, of which 3 GO terms were distributed in cellular 
component part (GO: 0005575, GO: 0043231 and GO: 
0043227). Furthermore, more highly ranked terms for the targets 
genes with higher expression were both in biological process, 
based on the significant analysis of GO terms between moso 
bamboo and ma bamboo. For example, on account of the 
category of GO terms, the ratio of biological process, cellular 
component and molecular function in moso bamboo was 40%, 
30% and 30% of the total of GO terms, respectively, while the 
corresponded ratio were change to 65.6%, 12.5% and 21.9% in 
ma bamboo. Higher proportion of biological process in ma 
bamboo possibly reflected the significance of target genes, which 
might play more essential role than those in moso bamboo. 
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Validation of specific conserved miRNAs by qRT-PCR and 
sequencing 

Specific primers designed for mature miRNA were used to 
validate the expression of miR396, miR397, miR1432, and 
miR7748 in leaves of both moso bamboo and ma bamboo. The 
result of qRT-PCR demonstrated that miR396, miR397, 
miR1432, and miR7748 were expressed differendy in leaf of 
moso bamboo, among which miR396 was the highest one, 
followed by miR1432 and miR397, and miR7748 was the lowest 
one (Figure 5). 

While that in ma bamboo, miR397 was the highest one, 
followed by miR1432 and miR396, miR7748 was not detectable 
(File S8). To prove the existence of miRNAs in ma bamboo, stem- 
loop primers were also directly used to amplify the mature miRNA 
sequences using specific templates from moso bamboo and ma 
bamboo, respectively, and the amplified products were sequenced. 
The result showed the sequences of mature miRNA in moso 
bamboo were all consistent with the sequences of miR396, 
miR397, miR1432, and miR7748 correspondingly. However, only 
three pairs of primers had specific amplification products for ma 
bamboo, one is the sequence of miR396,consistented with that of 
moso bamboo, while the other sequences were mismatched with 
miR397 and miR1432, respectively (File S9), which indicated the 
qRT-PCR result of miR397 and miR1432 for ma bamboo is false 
positive. 

Moreover, the primers designed for precursors of miR396, 
miR397, miR1432, and miR7748 were used to amplify with DNA 
templates of moso bamboo and ma bamboo, respectively. The 
sequencing result demonstrated the sequences of moso bamboo 
were matched with the genomic sequences and could form the 
stem-loop structures (File S10). However, there were no specific 
PCR products from ma bamboo. In order to identify the targets of 
these miRNAs, we had screened the published ma bamboo data of 
transcriptome [42] and ESTs [43] in NCBI and bambooGDB 
database [44] using the targets of miR396, miR397, miR1432, 
and miR7748, respectively, and no matched sequences was found 
for them, which indicated that the target genes in moso bamboo 
and ma bamboo are not comparable. On the basis of these results, 
miR397, miR1432, and miR7748 were specific conserved in the 
leaf sample of moso bamboo. 

Conclusions 

To the best of our knowledge, this work presents the firsdy 
comparative profiling of microRNAs in representative monopodial 
bamboo (moso bamboo) and sympodial bamboo (ma bamboo), 
based on high-throughput sequencing analysis of small RNA. A 
large number of conserved miRNAs, novel miRNAs and target 
genes were identified, which are worthy of further investigation to 
improve our understanding of the regulatory mechanisms of 
miRNA in biological processes in bamboo. Target prediction of 
differential miRNAs indicated that SPL genes are mainly regulated 
by miRNAs. Functional annotation of genes targeted by differen- 
tial miRNAs revealed that high significant biological processes are 
transcription regulation and photosynthesis. The expression of 
novel miRNAs and target genes involved in developmental stage 
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